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Understanding the damage of DNA bases from hydrogen abstraction by free OH radicals 
is of particular importance to reveal the effect of hydroxyl radicals produced by the sec- 
ondary effect of radiation. Previous studies address the problem with truncated DNA bases 
as ab-initio quantum simulation required to study such electronic spin dependent processes 
are computationally expensive. Here, for the first time, we employ a multiscale and hybrid 
Quantum-Mechanical-Molecular-Mechanical simulation to study the interaction of OH radicals 
with guanine-deoxyribose-phosphate DNA molecular unit in the presence of water where all the 
water molecules and the deoxyribose-phosphate fragment are treated with the simplistic classical 
Molecular-Mechanical scheme. Our result illustrates that the presence of water strongly alters the 
hydrogen-abstraction reaction as the hydrogen bonding of OH radicals with water restricts the 
relative orientation of the OH-radicals with respective to the the DNA base (here guanine). This 
results in an angular anisotropy in the chemical pathway and a lower efficiency in the hydrogen 
abstraction mechanisms than previously anticipated for identical system in vacuum. The method 
can easily be extended to single and double stranded DNA without any appreciable computational 
cost as these molecular units can be treated in the classical subsystem as has been demonstrated here. 
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PACS numbers: 

It is known that the megavoltage radiation (X/7-rays, 
a-particles, and heavy ions) ionize the water molecule 
and create neutral OH-radicals [1^ . Various effects of the 
ionizing radiation on the biological systems that ranges 
from the development of genetic aberrations, carcinogen- 
esis to aging have attributed to the role of free radi- 
cals |2J. In particular OH-radicals are major contribu- 
tors in the single/double strand breaking of the DNA 
molecules as the 2/3 of the surrounding environment of 
the DNA molecules in the cell-nucleus is composed of wa- 
ter molecules. OH-radical, a polar diatomic molecule is 
highly reactive. They are electrically neutral with a mag- 
netic moment that is produced by nine electrons with an 
unpaired spin in the outermost open shell. In-vivo, they 
have a very short life-time that is reported within a few 
nano-seconds [3 . Because of their short life-time, OH 
radicals generated within 1 nm from the surface of the 
DNA molecule can remove a hydrogen ion and form a 



water molecule. The high reactivity of radical OH is at- 
tributed to the possible pairing of opposite spin electrons 
from OH and the target molecule. Due to reactivity to 
solvent accessibility, a free radical reaches a DNA active 
site and removes a hydrogen ion from sugar moiety (e.g., 
a deoxyribose) [IHS] or nucleobase . 

The motion of the electrically neutral OH-radical is 
governed by the thermal diffusion process. Like any di- 
atomic molecule, their energy can be expressed as the 
sum of electronic, vibrational, translational and rota- 
tional parts [9 . The conformation of the molecular levels 
with the environment during the dynamical motion of the 
molecule allows the molecule to find its optimal chemical 
pathway. From the point of view of the tetrahedral con- 
formation of the electronic charge distribution in Oxygen, 
the hydrogen abstraction pathways by OH are expected 
to show angular anisotropy with a peak around 109.5 
degrees for H-OH angle. Recent simulations of the hy- 
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drogen abstraction in the vacuum by OH-radical pTHTS] 
also observe that the occurrence of the chemical reac- 
tion strongly depends on the relative orientation of the 
OH radical and the H atom of the donor molecule. In 
vacuum simulation, the radical OH is quasi free to ori- 
ent itself appropriately with guanine hydrogens. How- 
ever, in solution, the OH radicals are expected to form 
hydrogen bonding with water molecules which would in 
turn restrict their relative orientations with respect to 
guanine-hydrogens and we may see a retardation in hy- 
drogen abstraction. Thus, understanding the role of wa- 
ter molecules for quantitative analysis of the DNA initial 
damage by ionizing radiation in a realistic environment 
is necessary. However, the exponential growth of com- 
putational procedure with the system size is the main 
challenge in performing realistic calculations. A realistic 
environment contains of thousands of water molecules 
and is too expensive to be handled by any ab-initio cal- 
culation - which is required to describe the H-abstraction 
reaction. The multiscale, hybrid Quantum-Mechanical- 
Molecular-Mechanical (QM/MM)) method is a practical 
approach that allows significantly lower the computa- 
tional effort, yet suitable to address the electronic re- 
arrangement of the active site and consequent breaking 
and forming of chemical bonds. In a hybrid QM/MM, 
a small subsystem containing the active site of the re- 
acting molecules is treated with an ab-initio QM method 
while the rest of the system including the solvent, if any, 
is treated with a classical Molecular-Mechanical (MM) 
force-field. Application of a QM method is essential for 
the calculation of the redistribution of electronic charges 
of the active site on-the-fiy, as needed for the conforma- 
tions of the dynamical trajectory to simulate the chemical 
reaction pathways. 

In order to address the issue of hydrogen abstraction 
from a realistic system of a DNA guanine base [TOHIS], 
we consider a guanine base coupled with a deoxyribose 
and a phosphate backbone truncated at both ends of a 
single stranded DNA base unit - we refer it as a DNA 
residue which is one oxygen short (the 03' from the ad- 
jacent residue) to that of a nucleotide. During QM/MM, 
the deoxyribose and the phosphate group are treated in 
the MM subsystem and hence this truncation is not ex- 
pected to have any direct consequences in the electronic 
distribution of the system and hence in the hydrogen ab- 
straction process. However, it would certainly have trun- 
cation effect coming through the perturbation of the QM 
system by MM charges as we will always have with any 
truncated system. Our overall system consists of this 
DNA residue dissolved in water with one water molecule 
close to the guanine NH2 group converted in to a rad- 
ical OH by deleting a hydrogen. The solvated system 
was first given a short energy minimization using OPLS 
force field [24l [25] so as to eliminate any bad contacts 
with water molecules arising from solvation. The OH- 
radical seems to form hydrogen bonding with the NHi 



and NH2 groups of guanine. To understand the hy- 
drogen abstraction mechanism, the system is then sub- 
jected to the hybrid GROMACS-CPMD [22] Quantum- 
Mechanical-Molecular-Mechanical (QM/MM) simulation 
protocol where the active site consisting of the OH and 
the guanine base are treated in the QM subsystem and 
the deoxyribose, the truncated phosphate backbone and 
water are treated in the MM subsystem. As the deoxyri- 
bose and the truncated phosphate backbone are treated 
in the classical MM model, and they lie in close viscinity 
of the active site treated in QM we do not use end capping 
atoms for the phosphate unit as they will bring additional 
point charges close to the QM subsystem. The QM sub- 
system is optimized with a plane-wave based DFT model 
as implemented in CPMD. The QM/MM partition is be- 
ing made at the Ng— C bond connecting the guanine ni- 
trogen to the deoxyribose. At the guanine-deoxyribose 
N-C bond-cut the link-atom prescription of QM/MM is 
used to saturate the QM subsystem for the QM simula- 
tion. The rest of the system - the deoxyribose, the phos- 
phate, and the water are treated in the MM subsystem 
using the Molecular Mechanics force field OPLS [24l [25] 
as implemented in GROMACS [26 . 

The overall computational box size is taken as 
(42.22A X 49.I4A X 40.78A) which contains the whole 
system: the guanine-based DNA residue (guanine- 
deoxyribose-phosphate unit), an OH radical, and 2811 
TIP3P [27 water molecules. A subset of the system 
containing the DNA residue and water within about 5A 
is shown in Figure la. For the QM subsystem, which 
contains the guanine and radical OH, a cubic cell of 
size (14.82A x 18.52A x 10.58A) is used together with 
Poisson solver of Martyna and Tuckerman [30 for the 
wave function optimization in CPMD. For both MM and 
QM/MM, we use a cutoff of lOA for the Coulomb and van 
der Waals interactions. In GROMACS-CPMD QM/MM, 
CPMD provides the optimized wave function and forces 
for the QM subsystem which is appropriately perturbed 
by the MM subsystem and GROMACS controls Molecu- 
lar Dynamics the simulation. For our QM/MM, we have 
performed a constant pressure and temperature (NPT) 
MD simulation with a reference temperature of 300K 
and a time step of Ifs. Our CPMD wave-function op- 
timization is implemented in a plane-wave basis within 
local spin density approximation (LSDA) with an en- 
ergy cutoff of 40 Rydberg (Ry), and with Becke \2E\ 
exchange and Lee- Yang-Parr (BLYP) gradient-corrected 
functional [29 . Norm conserving ultrasoft Vanderbilt 
pseudo-potentials were used for oxygen, hydrogen, nitro- 
gen and carbon. Three different sets of simulation were 
performed - two in solution and one in vacuum. For gua- 
nine spin singlet (5 = 0), which corresponds to a doublet 
multiplicity for the guanine- OH system, two simulations 
were performed - one with solvent water and the other 
without solvent. For guanine spin triplet {s = 1), which 
could lead to a quartet multiplicity for the guanine- OH 
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system, only one simulation has been performed with sol- 
vent. This is to be mentioned here that for this latter sys- 
tem, an earlier simulation in vacuum [I3l [14] performed 
using CPMD, found no hydrogen abstraction from gua- 
nine. The results presented in this work are representa- 
tive of several runs of QMMM molecular dynamics, start- 
ing from a large number of possible initial conditions. 

Our hybrid QM/MM simulations consist of two classes 
of spin-restricted calculations, as the total spin along the 
quantum axis is subjected to the constraints Sz = 1/2, 
and 3/2, corresponding to doublet and quartet spin con- 
figurations. In both spin-state calculations, same initial 
configurations were used for the guanine-OH subsystem 
which was taken from an OPLS force-field generated en- 
ergy minimized structure of the total system in the pres- 
ence of solvent water. To understand the effect of solvent 
on hydrogen abstraction, we perform QM/MM molecu- 
lar dynamic simulations of the MM energy minimized 
structure of the solvated system (guanine-deoxyribose- 
phosphate-OH) obtained with OPLS force-field. The 
initial configuration here corresponds to H-OH distance 
of 1.92 A and H-O-H angle of 108.92 degrees with the 
NH group hydrogen The initial MM energy minimized 
configuration is characterized by hydrogen bonding of 
OH with NH and NH2 groups of guanine and as well 
as with one solvent water. Where the H-bond lengths 
for Htvi— OH, Htvs— OH, and Hsoi—OH are given by 
1.92 A, I.97A, 1.68 A, respectively. The initial bond angle 
for Hat.-O-H is 108.42 degree and that of Hats-O-H is 
153.95 degrees. Consequently, the initial configuration, 
which is obtained by MM energy minimization, places 
Hati in a favorable position to form water with the rad- 
ical OH (since Hati— O— H angle is 108.42 degree which 
is close to the tetrahedral angle of 109.5 degrees). The 
same initial structure has been used to study the vacuum 
property by simply deleting the water from the solvated 
energy minimized configuration. Using a similar system 
in vacuum with a H-OH distance of 1.5 A, it has been 
reported earlier J3l [14] that the hydrogen abstraction 
does depend on the spin of guanine. For guanine spin 
singlet (guanine-OH doublet), the hydrogen abstraction 
is found to start in about 0.07 ps in vacuum and in 0.1 
ps in solution. For a time step of 2 fs, the simulation ex- 
hibits numerical instability with the abstracted hydrogen 
is found to have strong kinetic energies leading to bounc- 
ing back to guanine and then lost to water. As the time 
step is reduced to 1 fs, numerical stability returns back 
and we notice a gradual process of hydrogen abstraction 
from guanine. Even with 1 fs time step, the kinetic en- 
ergy of the abstracted hydrogen remains high. For simu- 
lation of the system in vacuum, this KE is found to get 
translated into the KE of the of the newly formed water 
and the resulting water molecule is found to move away 
quickly from the guanine. However, in the presence of sol- 
vent, the translational motion of OH is being restricted 
by its hydrogen bonding with solvent water molecules 
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FIG. 1: A smaller subset of the complete system containing a 
DNA residue with guanine base surrounded by about 5A layer 
of water is shown in la. In lb, the water molecules which 
are hydrogen bonded with the DNA residue and with the 
OH which are directly or indirectly influencing the hydrogen 
abstraction process are shown. Carbon, oxygen, nitrogen, 
phosphorus and hydrogen atoms are shown as green, red, blue, 
gold and white, respectively. 

and thus we see quite a lot of fluctuation in bond length 
and angle in the newly formed water for quite some time. 
For guanine triplet (guanine-OH quartet), we continued 
the QM/MM simulation for 1.5ps (about times the time 
frame needed to see hydrogen abstraction in singlet) but 
there was no sign of any significant H-OH attraction be- 
tween the radical OH and the amino group guanine hy- 
drogen as was expected from spin consideration. 

In Fig. [ija), we show the DNA residue dissolved in 
water and in Fig. [ijb) we show DNA residue with some 
of the explicit hydrogen bonding contacts with water 
molecule which directly or indirectly influence the hydro- 
gen abstraction dynamics. The QM/MM Molecular Dy- 
namics (MD) simulations have been performed for 1.5 ps 
and the MD trajectories are studied in VMD [31 . Car- 
bon, oxygen, nitrogen, phosphorus and hydrogen atoms 
are shown with colors green, red, blue, gold and white, 
respectively. According to our results, a dehydrogenation 
of the nucleotides takes place for a system with Sz — 1/2 
(total spin-doublet). Fig. |2] shows the time sequence of 
the hydrogen abstraction from Ni-site of guanine and the 
hydrogen bonds formation during the process that leads 
to the formation of a water molecule. 

Fig. [sja) shows the oxygen atom distances from the 
OH radical to the Hni in the spin singlet guanine in 
vacuum (filled red circles). As it is seen the abstrac- 
tion of Hydrogen occurs around t ^ 0.07 ps. In solution 
(open black circles), however, because of the hydrogen 
bond network forming with water molecules in the solu- 
tion and the OH-radical (see Fig.jl]) the time for hydrogen 
abstraction is longer and it occurs at t ~ 0.12 ps. This 
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FIG. 2: Time sequence of the Hydrogen abstraction at t = 0, 
t = 0.1 ps, and t = 0.5 ps are shown for simulations in vacuum 
(a,b,c) and in the presence of water (d,e,f). The Hydrogen 
bonds formed between water molecules in solution and the 
OH-guanine molecules is shown by dashed lines. The numbers 
in the figures denote the bond distances. The hydrogen bonds 
elongate the time for hydrogen abstraction. 

is the main finding of this study as it reveals that the 
reaction in solution takes place with slower rate, approx- 
imately by a factor of two, and because of the ns lifetime 
of OH-radicals, roughly half of them have chance of mak- 
ing reaction. For comparison the spin triplet in solution 
(blue stars) are shown. A drop in i^7vi-0H seen in spin 
singlet state that is the indication of guanine dehydro- 
genation, disappears in case of spin triplet state where 
the spin-blockade effect is strong enough that it repels 
OH free radical. Clearly the hydrogen abstraction does 
not occur for triplet guanine in solution for the time scale 
of the present MD simulation. This confirms the previ- 
ous calculation that the injection of 5* = 1 exciton blocks 
the hydrogen abstraction [13l[14]. Fig.jsjb) and (c) show 
the evolution of the angle between Hni and OH radi- 
cal in the spin singlet guanine in vacuum (filled red cir- 
cles) and in solution (open black circles). The horizontal 
dashed line shows the angle 109.5 degrees. The vertical 
lines show the expected position of Hydrogen abstrac- 
tion, consistent with the distance-time calculation shown 
in Fig. [3]^ a). This confirms our hypothesis on the an- 
gular dependence of the hydrogen-abstraction of DNA 
molecule by free radicals, stated based on the circular 
symmetry of the diatomic molecule. Accordingly, the 
range of angles for which Hydrogen abstraction occurs is 
expected to be uniformly distributed in < (/? < 27r, but 
sharply distributed around 6 = 109.5, the angle between 
sp3 orbitals in oxygen prior to the water- formation. Af- 
ter hydrogen-abstraction, as seen in Fig.jsjb) and (c), the 
mean angle between two hydrogens in H2O molecule be- 
comes — 104.45. Small vibrations around = 104.45 is 
clearly visible. These vibrations are stronger in solution 



than in vacuum. Here (p and represent the spherical an- 
gles in a system of coordinates that the dipole-moment 
of OH-radical is along z-axis with oxygen at the center 
of coordinate. This is the reflection of the circular sym- 
metry of the diatomic molecule that strongly alters in 
aqueous solutions because of the hydrogen bonds between 
OH and neighboring water molecules. Figs. |3|b) and (c) 
clearly illustrate that the alignment of the diatomic OH- 
radical necessary for the hydrogen-abstraction is more 
difficult in solution because of the thermal fluctuations 
of the solvent that transfers to OH-radicals through the 
hydrogen-bonds formation with water molecules. 

Fig. [4] shows the Kohn-Sham energies (equivalent to 
potential energy in classical molecular dynamics) of the 
spin singlet and spin triplet of the guanine in the presence 
of the OH free radical as a function of time, calculated 
by the CPMD at T=300K corresponding to a canonical 
dynamics (constant temperature ensemble). A drop in 
Kohn-Sham energy, A^ks^ iii spin singlet multiplicity 
is indication of dehydrogenation of Uni in the guanine 
by OH free radical. Our calculation shows that AEks = 
-1.088 (-0.844) eV per OH-radical in solution (vacuum). 
Because |A^ks| (in — solution) > |A£^ks| (in — vacuum), 
the Hydrohen-abstraction in solution is energitically fa- 
vorable, however, because of the thermal fluctuations 
of the water molecules that disturb the motion of OH- 
radical, the reaction takes place with slower rate (note 
that the minimum of A£^ks is seen at a longer time in 
solution than in vacuum). Hence the network of Hy- 
drogen bonds between water molecules and OH-radical, 
studied in this work, is the bottle-neck of the Hydrogen- 
abstraction in solution. Finally, in the quartet spin con- 
figuration the repulsive exchange interaction, blocks the 
exchange of the hydrogen and hence the chemical reac- 
tion. 

In Figs.[5]the time evolution of c^h-oh as a function of 
initial distance between OH and guanine in solution and 
in vacuum are shown. The initial OH distance varies 
along a line to keep OH-Nhi and 011-Nh2 equidistance. 
This is a line perpendicular to a triangle edge that con- 
nects Nhi and Nh2 as shown in Fig. [ijb). The initial 
distances in Figs. [H] are (in-OH = 1-9, 2.3, 2.8A from 
(a) to (c). Comparing Figs, [sja) and (b), we observe 
that the increase of initial distance from (in-OH = 1.9 to 
2. 3 A does not alter the final product. Moreover the reac- 
tion time in both solution and vacuum appears not to be 
sensitive to the change of (in- oh within this range. How- 
ever, increasing initial distance to (in- oh = 2.8A lead to 
a dramatically different final state in solution where the 
abstraction of Nh2 is seen. This is mainly due to ef- 
fect of solution such that OH-radical gains translational 
and rotational energy from water molecules and crosses 
through transition state of 0H-A^//2- 

In summary, we have employed a multiscale Quantum- 
Mechanical-Molecular-Mechanical (QM/MM) simulation 
based on GROMACS-CPMD QM/MM protocol to study 
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FIG. 3: In (a) the evolution of the distances (in A) from 
oxygen atom in the OH radical to the Hni in the spin sin- 
glet guanine in vacuum (filled red circles), in solution (open 
black circles), and in the spin triplet in solution (blue stars) 
are shown. The hydrogen abstraction for singlet guanine in 
vacuum and solution occurs around t ^ 0.07 and t ^ 0.12 ps, 
respectively, and it does not occur for triplet guanine in solu- 
tion. In (b) and (c) the evolution of the angle between Hni 
and OH radical in the spin singlet guanine in vacuum (filled 
red circles) and in solution (open black circles) are shown. 
The horizontal dashed line shows the angle 109.5 degrees. The 
vertical lines show the expected position of Hydrogen abstrac- 
tion, consistent with the distance-time calculation shown in 
(a). 



the hydrogen abstraction by OH-radicals from a guanine- 
based DNA residue solvated in explicit TIP3P water. We 
found that the presence of water molecules slows down 
the process of hydrogen abstraction by OH radicals as 
the latter's thermal fluctuation is restricted by its hy- 
drogen bonding with solvent water molecules. To under- 
stand the effect of arbitrary separation of the OH rad- 
ical from target hydrogen on hydrogen abstraction, we 
repeated the simulation on two other H-OH distances. 
For arbitrary separations, the NH2 hydrogens also act as 
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FIG. 4: The Kohn-Sham energy of the QM-part as a func- 
tion of time and number of water molecules of spin-triplet 
DNA (top) and spin-singlet DNA (bottom) in the presence 
of one OH radical. A drop in Kohn-Sham energy is indica- 
tion of dehydrogenation of Ha^i in the guanine by OH free 
radical clearly seen in the spin-singlet DNA-OH system (bot- 
tom). The hydrogen abstraction in vacuum occurs around 
t ^ 0.1 ps. In solution this occurs around t ^ 0.2 ps, ap- 
proximately two times slower. The Hydrogen-bonds coupling 
between water molecules and DNA residue are responsible to 
strong fiuctuations in Kohn-Sham energy. 



potential targets for hydrogen abstraction as the ther- 
mal diffusion of OH plays a significant role in the hydro- 
gen abstraction dynamics of OH. The present simulation, 
performed on a larger system of guanine attached to a 
deoxyribose and phosphate, exhibits that hydrogen ab- 
straction occurs only when the guanine-H system is in its 
lowest spin doublet state; hydrogen abstraction is blocked 
for guanine- OH spin quartet. A similar effect of spin 
coupling on hydrogen abstraction obtained earlier with 
only guanine in vacuum. Present results, while comple- 
ment previous calculations on smaller system in vacuum, 
demonstrate that such multiscale simulation, employed 
for the first time to address hydrogen abstraction from 
DNA residue, could provide a basis to address the big- 



6 




0.03 



(b) 



0.06 0.09 

t(ps) 



0.12 0.15 



•-• in solution 
^ in vacuum 




0.03 0.06 0.09 0.12 0.15 

t(ps) 




0.06 0.09 

t(ps) 



0.15 



FIG. 5: The time evolution of c/h-oh as a function of initial 
distances du-ou = 1-9, 2.3, 2.8A. In (a) and (b) du-ou rep- 
resents the the distance between OH and Nhi as the hydrogen 
abstraction occurs from Nhi- In (c) hydrogen abstraction oc- 
curs from Nhi and Nh2 in vacuum and solution respectively. 



ger question of DNA damage by radicals in presence of 
solution. 
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